###################################################
## Demographic Analysis
## SCRIPT 4: Script to run mergesplit in the background
## at various VRA constraint strengths.
##
## 5/17/21
###################################################

library(ppmf)
library(geomander)
library(tidyverse)
library(sf)
library(redist)

la_all <- readRDS("../../data/LA/la.Rds") %>%
  rename(sld_upper = sld_up,
         sld_lower = sld_low) %>%
  filter(sld_lower != 106)

la_map <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = pop,
    pop_tol = 0.05,
    existing_plan = sld_lower
  )

la_map4 <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = v4_pop,
    pop_tol = 0.05,
    existing_plan = sld_lower
  )

la_map12 <- la_all %>%
  st_set_crs(26915) %>% 
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = v12_pop,
    pop_tol = 0.05,
    existing_plan = sld_lower
  )

la_map19 <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = v19_pop,
    pop_tol = 0.05,
    existing_plan = sld_lower
  )


##################################################################
## useful constant objects
district_focus <- c(29, 61, 63:71, 101, 18, 62, 60)

ebr_map <- la_map %>% filter(sld_lower %in% district_focus)
ebr_map4 <- la_map4 %>% filter(sld_lower %in% district_focus)
ebr_map12 <- la_map12 %>% filter(sld_lower %in% district_focus)
ebr_map19 <- la_map19 %>% filter(sld_lower %in% district_focus)
##################################################################
## plans

strength_vals <- c(0, 1, 10, 100, 500, 1000, 10000, 25000, 50000)

for (i in 1:length(strength_vals)) {
  
  plans_smc <- redist_mergesplit(ebr_map,
                                 init_plan = match(ebr_map$sld_lower,
                                                   unique(ebr_map$sld_lower)),
                                 nsims = 50000,
                                 constraints = list(vra = list(strength = strength_vals[i],
                                                               min_pop = ebr_map$pop_black)),
                                 warmup = 1) %>%
    add_reference(ebr_map$sld_lower, "2010 Adopted")

  filename <- paste("../../data/LA/sim/vra_test/ms_", strength_vals[i], ".rds", sep = "")
  saveRDS(plans_smc, file = filename); rm(plans_smc)

  plans4_smc <- redist_mergesplit(ebr_map4,
                                  init_plan = match(ebr_map4$sld_lower,
                                                    unique(ebr_map4$sld_lower)),
                                  nsims = 50000,
                                  constraints = list(vra = list(strength = strength_vals[i],
                                                                min_pop = ebr_map$v4_pop_black)),
                                  warmup = 1) %>%
    add_reference(ebr_map4$sld_lower, "2010 Adopted")

  filename <- paste("../../data/LA/sim/vra_test/ms4_", strength_vals[i], ".rds", sep = "")
  saveRDS(plans4_smc, file = filename); rm(plans4_smc)

  plans12_smc <- redist_mergesplit(ebr_map12,
                                  init_plan = match(ebr_map12$sld_lower,
                                                    unique(ebr_map12$sld_lower)),
                                  nsims = 50000,
                                  constraints = list(vra = list(strength = strength_vals[i],
                                                                min_pop = ebr_map$v12_pop_black)),
                                  warmup = 1) %>%
    add_reference(ebr_map12$sld_lower, "2010 Adopted")

  filename <- paste("../../data/LA/sim/vra_test/ms12_", strength_vals[i], ".rds", sep = "")
  saveRDS(plans12_smc, file = filename); rm(plans12_smc)
  
  plans19_smc <- redist_mergesplit(ebr_map19,
                                   init_plan = match(ebr_map19$sld_lower,
                                                     unique(ebr_map19$sld_lower)),
                                   nsims = 50000,
                                   constraints = list(vra = list(strength = strength_vals[i],
                                                                 min_pop = ebr_map$v19_pop_black)),
                                   warmup = 1) %>%
    add_reference(ebr_map19$sld_lower, "2010 Adopted")
  
  filename <- paste("../../data/LA/sim/vra_test/ms19_", strength_vals[i], ".rds", sep = "")
  saveRDS(plans19_smc, file = filename); rm(plans19_smc)
}
